SMALLFFT.EXE   eine Abgemagerte Version um die Programmierung  
SMALLFFT.C     besser verstehen zu k”nnen!

FFT.ZIP        die FFT-Routinen, alt (1988) aber gut!
	       jedoch mangelhaft dokumentiert...




	       Kurzbeschreibung FFT.ZIP
---------------------------------------------------------------------------
die aufgenommenen Werte werden im real[] Array abgelegt,
das imag[] Array wird mit Nullen aufgef�llt...

dann werden die Globalen Variablen samples und power gesetzt:
samples muá ein Potenz von 2 sein: 2,4,8,16,32,64,128,256,512,1024,2048,...

int     samples, power;
double  real[2048], imag[2048], max;
samples = Anzahl aufgenommener Werte im real[]/imag[] Array;
power =  log10((double)samples) / log10((double)2.0);

Jetzt wird die Funktion fft() ausgef�hrt!
dieses ist sehr schnell, extrem viel schneller als die normale FT.
das Real/Imagin„r-Spektrum befindet sich danach im real[]/imag[] Array.
Allerdings recht merkw�rdig durcheinandergew�rfelt (das ist bei FFT normal)!
die permute(n) liefert uns die Position des Spektrums n in real[]/imag[]
Beispiel: 
		      t=permute(n) 
		      r=real[t]
		      i=imag[t]
		      Betrag=sqrt(r*r + i*i) 
		      return Betrag;

dieses Beispiel ist die Funktion magnitude(n), welche den Betrag liefert!
Wenn wir uns eh nur f�r das Betragspektrum interessieren benutzen wir
also magnitude() und haben mit permute() nix mehr am Hut!



Iss auch ne R�ckw„rts (inverse) Fast-Fourier-Transfo m”glich?
Aber klaro, laut dem schlauen FFT-Transformation Buch das ich mal kurz
in den Fingern hatte ist dazu lediglich der Imagin„rteil zu invertieren, 
(also mit -1 zu multiplizieren), fft() aufrufen, und dann wieder den
Imagin„rteil invertieren.
aber dabei ist die Richtige Reihenfolge im Array wichtig (siehe permute)!

fft();   
for(i=0;i<samples;i++) imag[i]*=(-1);
fft();
for(i=0;i<samples;i++) imag[i]*=(-1); 

nach diesem kurzen Programm stehen nicht etwa wieder die urspr�nglichen
Werte im Array, sondern Datenm�ll, weil die Verw�rfelung ignoriert wurde!
Richtiger w„re:

fft();   
for(i=0;i<samples;i++)
		      {
		       t=permute(i)
		       real2[i]=real[t];
		       imag2[i]=imag[t];
		      }
for(i=0;i<samples;i++) real[i]=real2[i], imag[i]=imag2[i];
for(i=0;i<samples) imag[i]*=(-1);
fft();
for(i=0;i<samples) imag[i]*=(-1); 
for(i=0;i<samples;i++)
		      {
		       t=permute(i)
		       real2[i]=real[t];
		       imag2[i]=imag[t];
		      }

Nach der Ausf�hrung dieses Programs stehen die ursp�nglichen Werte mit
kleinen/grossen Fehlern wieder im real[]/imag[] Array und im 
Zwischenspeicher-Array real2[]/imag2[] finden wir das Real/Imagin„rspektrum!


welches Format hat denn nu dasch Real/Imagin„r/Betragsspektrum?
0 = Gleichanteil ....... samples-1 = samplerate
leider k”nnen wir das Array nur bis zur H„lfte (samples/2 - 1) nutzen weil
ein gewisser Nyquist das sogenannte Nyquist-Kriterium aufgestellt hat,
nach dem bei Abtastung maximal Frequenzen der halben Abtastrate abgetastet
werden k”nnen. 
Erkl„rung mittels Systemtheorie: (bahh, w�rg, igittigitt)
Ich versuche es unmathematisch visual anschaulich zu erkl„ren,
also Zeichblatt und Stift rauskramen und versuchen das Folgende zu zeichnen:

wir betrachten die Abtastung im Freqenzbereich:
die kontinuierliche Abtastung (Kammfunktion im Frequenz- und Zeitbereich)
mit Abstand   samplerate (Hz) im Freqenzbereich und
mit Abstand 1/samplerate (s)  im Zeitbereich
wird im 
Freqenzbereich gefaltet mit dem Eingangssignalspektrum und im
Zeitbereich multipliziert mit der Eingangsignal.

Faltung ist normalerweise besch... zu zeichnen aber wegen den DiracstӇen
und der Kammfunktion wirds viel einfacher, wir brauchen bloss das
Eingangssignalspektrum um jeden Diracstoá zu zeichnen, als wenn jeder
Diracstoá f=0Hz h„tte, wenn wir das gezeichnet haben, 
verstehen wir nicht nur das Nyquist-Kriterium, 
wir wissen auch wie unser Spektrum zwischen samples/2 ... samples aussieht!
Das sind n„hmlich nicht weiter als die negativen Frequenzen des um
Samplerate gefalteten Eingangssignalspektrums.

Also entspricht der linke Rand (0) und der rechte Rand (samples)  des Arrays 
(permute beachten) dem Gleichanteil (0 Hz) des Eingangsspektrums, da der
maximale adressierbare Wert im Array jedoch samples-1 ist finden wir den 
Gleichanteil nur bei real[permute(0)],imag[permute(0)]... 


Realteil:
das Signal von 0...samples/2 - 1 wird um samples/2 gespiegelt.
Imagin„rteil:
das Signal von 0...samples/2 - 1 wird um samples/2 gespiegelt, und zus„tzlich
invertiert (also das Vorzeichen gewechselt).

Betragspektrum: um das Betragsspektrum in den Zeitbereich zu transformieren
zerlegen wir es in eine gerade und eine ungerade Funktion,
der gerade Anteil kommt nach real[]
der ungerade Anteil kommt nach imag[]
dann wird Realteil und Imagi„rteil wie angegeben gespiegelt!

f�r die inverse FFT ist nat�rlich noch die weiter oben erw„hnte doppelte
Invertierung des Imagin„rteils notwendig!!!
-----------------------------------------------------------------------------


Informationstragend ist dabei eigendlich nur der Bereich von 0...samples/2-1
entsprechend f=0 bis f=samplerate/2 also von 0 Hz bis zur halben Abtastfrequenz.

Wie meine Experimente mittels Soundkarte und Funktionsgenrator zeigten,
ist beim Realspektrum und beim Imagin„rspektrum das Vorzeichen dauernd am
kippen, schade...
das Betragsspektrum wird mittels Arctan(imag[permute(i)]/real[permute(i)]) 
berechnet, leider eine ungerade und damit vorzeichensensitive Funktion...

das Betragsspektrum ist nicht nur sehr aussagekr„ftig, es ist auch wegen
der Quadratbildung nicht vorzeichensensitiv also stabil!!!

es gibt da allerdings noch ein paar Probleme auf die ich gestoáen bin:
bedingt durch die Abtastung kommt es wenn sich ein Vielfaches der 
Abgetasteten Frequenz sich der Abtastrate n„hert zu einer Schwebung
(Amplitudenmodulation), daraus resultieren mit Sicherheit Fehler!

die FFT hat auch noch ein „hnlich gelagertes Problem:
wie wir aus unserer Zeichnung ersehen k”nnen, ist die diskrete 
Fourier Transformation bedingt durch die Abtastung 
im Freqenzbereich periodisch, daraus ergibt sich der Umkehrschluá das
auch irgendeine Art von Periodischen Verhalten im Zeitbereich vorliegt!
leider trifft dieser Umkehrschluá zu !!!!!!!
wenn irgendein Vielfaches der Periodenl„nge der Abgetasteten Funktion sich
der Zeit (samples/samplerate) also (Anzahl_Werte/Abtastfrequenz) n„hert
dann werden die Spektralinien hoch und schmal, wenn dies nicht zutrifft
werden die Spektallinien entsprechend der Nicht�bereinstimmung flach und
Breit, die Energie liegt also in der Fl„che!

Auáerdem steht zu bef�rchten das sich die Beiden Probleme �berlagern
(belagern) weil sie n„hmlich unter den selben Bedingungen auftreten.
M”glicherweise handelt es sich hierbei nicht um eine šberlagerung zweier
Probleme sondern lediglich um ein von Zeitbereich in den Frequenzbereich 
transformiertes Problem, verwandt sind bei Probleme auf jeden Fall!




Es gibt zwei m”glich L”sungen dieses Žrgerlichen Problems,

1. man benutzt eine geeignete Fensterfunktion windows(), daá heiát:  
   man multipliziert das Signal im Zeitbereich mit einer sanft Anklingenden
   und sanft ausklingenden Funktion (sin,gauss...) so das zum Beispiel 
   eine Abgetastete Funktion nicht schlagartig anf„ngt und abbricht
   sondern langsam ein und wieder ausgeblendet wird.

2. man nimmt das abzutastende Signal mehrfach mit unterschiedlichen 
   Abtastraten auf und �berlagert die entstehenden Betragsspektren!
   dabei empfehle ich instinktiv das Quadratische Mittel zur šberlagerung
   heranzuziehen. Nat�rlich m�ssen die Frequenzen umgerechnet werden,
   und genau da liegt ein weiterer gef„hrlicher Fallstrick, denn
   die Mega-Blaster 16 (Mozart) Soundkarte emuliert ja nur eine 
   Soundblaster Pro Soundkarte, das heisst das sich aufgrund der 
   unterschiedlichen Hardware und des Begrenzten Wertebereich (8.Bit) 
   der Timerconstante kleine Unterschiede betreffend die tats„chliche 
   Abtastfreqenz ergeben. Die Folgen f�r ein šberlagerung w„ren Katastrophal,
   wir m�ssen uns also zus„tzlich die M�he machen die tats„chlichen
   Abtastfrequenzen mittels der weiter oben erw„hnten Schwebung zu ermitteln!
-----------------------------------------------------------------------------


ach ja, mittels der Funktion sbrec(samples); und sbrec2(samples); werden
Samples Bytes aufgenommen, jeweils im Normal- b.z.w. Highspeed-Mode.

if(HighSpeed) sbrec2(samples); else sbrec(samples);

der Highspeedmode sollte keinesfalls auf Frequenzen <13kHz angewendet werden!
der Normalmode kann maximal bis 15kHz/22kHz Samplerate verwendet werden.
die Mega-Blaster 16 (Mozart) hat eine maximale Abtastrate vom 48KHz.
praktisch k”nnen mit dieser Soundkarte Signale von 10Hz bis 21kHz aufgenommen
werden, ein linearer Freqenzgang ergibt sich von 50Hz bis 18kHz...



-----------------------------------------------------------------------------
Die VGA-Grafik Programmierung erfolgt aus Geschwindigkeitsgr�nden direkt, 
also nicht �ber die Borland_C Grafikbibliothek.
Um meine schwerverdauliche Highspeed-Grafik besser (�berhaupt) zu verstehen 
empfehle ich dringend(zwingend) die Lekt�re des folgenden Buches auf
dem meine Grafik-Routinen beruhen, es handelt sich allerdings nicht bloá um 
eine Abschrift der dort aufgef�hreten Routinen (die sind n„hmlich in Pascal)
sondern um auf den Informationen dieses Buches basierenden Code !!!

	 -------------------------------------------------
	 |  Die Programmierung der EGA/VGA Grafikkarte   |
	 |  Matthias Uphoff                              |
	 |  TWJ 132                                      |
	 |  ISBN 3-89319-274-3                           |
	 -------------------------------------------------

Das Buch ist in deutscher Sprache geschrieben,und steht in unserer Bibliothek.
Die Diskette braucht ihr euch gar nicht erst ausleihen.    (�berf�ssig!)